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Abstract 

We show that a simple and straightforward rational approximation to the Thomas- 
Fermi equation provides the slope at origin with unprecedented accuracy and that 
relatively small Pade approximants are far more accurate than more elaborate ap- 
proaches proposed recently by other authors. We consider both the Thomas-Fermi 
equation for isolated atoms and for atoms in strong magnetic fields. 



1 Introduction 

The Thomas-Fermi (TF) equation has proved useful for the treatment of many 
physical phenomena that include atoms [1-5], molecules [3,6], atoms in strong 
magnetic fields [1,4,5], crystals [7] and dense plasmas [8] among others. It is 
well-known that an accurate solution to that equation is based on the accurate 
calculation of the slope at origin [9-11]. In particular we mention the rational 
approximation in terms of Pade approximants [12,13] because it is relevant to 
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present discussion. As expected from the great physical significance of the TF 
equation, its mathematical aspects have been studied in detail [14-17]. 

In spite of being a quite old problem in theoretical physics, there has re- 
cently been a renewed interest in analytical solutions to the TF equation. For 
example, Liao [18] and later Khan and Xu [19] proposed the application of 
the so called homotopy analysis method (HAM). More recently Parand and 
Shahini [20] showed that a pseudospectral method based on Chebyshev poly- 
nomials is more accurate than HAM. Unfortunately, Parand and Shahini [20] 
were not aware that we had earlier shown that the well-known straightforward 
Pade approximants are much more accurate than HAM [21]. Therefore, they 
compared their pseudospectral results with the rather insufficiently accurate 
HAM ones that do not constitute a suitable benchmark. 



The purpose of this paper is twofold: first we compare the results of the simple 
and straightforward Hankel-Pade method (HPM) [21] with the supposedly 
accurate Chebyshev pseudospectral ones [20]. Second, we show that the HPM 
also gives accurate results for the TF equation for atoms in strong magnetic 
fields that has not been treated before in this way. 

In Section [2] we briefly introduce the TF equations for both problems and show 
how to transform them into more tractable differential equations. In Section [3] 
we outline the main ideas behind the HPM, apply it to both TF equations, 
and compare HPM and Chebyshev pseudospectral results for isolated atoms. 
In Section H] we summarize the main results and comment on other approaches 
to nonlinear ordinary differential equations. 
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2 The Thomas— Fermi equation 



The dimensionless form of the TF equation for atoms [1-5] 



'u(x) 3 



u"(x) = \ - LJ -, u(0) = 1, u(oo) = 
V x 



is an example of two-point nonlinear boundary- value problem. When solving 
this ordinary differential equation one faces the calculation of the slope at 
origin u'(0) that is consistent with the physical boundary conditions indicated 
in equation (CQ). 

In order to make this letter self contained we review the main features of the 
HPM [21]. It is convenient to define the variables t = x 2 and f(t) = -u(i 2 ) 1//2 , 
so that the TF equation becomes 

f(t)f(t) + f\tf] - f(t)f(t) - 2t 2 /(t) 3 = (2) 
We expand the solution f(t) in a Taylor series about t = 0: 

oo 

/(*) = E/i** ( 3 ) 

3=0 

where the coefficients fj depend on the only unknown one ji = /"(0)/2 = 
u'(0)/2. On substituting the series into equation (j2J) we easily calculate as 
many coefficients fj as desired; for example, the first ones are 

fo = 1, fl — 0, /3 = -, f± = — , /5 = — , ••• (4) 

The TF equation for atoms in strong magnetic fields is somewhat simpler 
[1,4,5] 

u"{x) = \Jxu(x), w(0) = 1, u{oo) = (5) 
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By means of the same change of variables we obtain 

t [f(t)f"(t) + f{tf] - f(t)f'(t) - 2t 4 f(t) = (6) 

as well as the solution in the form of the Taylor series ([3]) with coefficients 
/o = l, /i = 0, / 3 = 0, /4 = -y, / 5 = ^,... (7) 

where fi = u'(0)/2 as in the preceding case. 
3 The Hankel-Pade method 

The HPM is based on the transformation of the power series ([3]) into a rational 
function or Pade approximant 

\MIN](t) = (8) 

One would expect that M < N in order to have the correct limit at infinity; 
however, in order to obtain an accurate value of it is more convenient to 
choose M = N + d, d = 0, 1, . . . as in previous applications of the approach 
to the Schrodinger equation (in this case it was called Riccati-Pade method 
(RPM)) [22-30]. 

The rational function (jS]) has 2N + d + 1 coefficients that we may choose so 
that T([M/N],t) = 0(t 2N+d+1 ), where T(f,t) = stands for either equation 
(j2J) or (jUJ), and in both cases the coefficient ji remains undetermined. If we 
require that T([M/N],t) = 0(t 2N+d+2 ) we have another equation from which 
we obtain f 2 . However, it is simpler and more practical to proceed in a different 
(and entirely equivalent) way and require that 

[M/N](t)~ f 3 t j = 0(t 2N+d+2 ) (9) 

3=0 
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In order to satisfy this condition it is necessary that the Hankel determinant 
vanishes 



Hf) - I 1^=0,1,.. JV ~~ °' ( 10 ) 

where D = N + 1 is the dimension of the Hankel matrix. Each Hankel deter- 
minant is a polynomial function of f 2 and we expect that there is a sequence 
of roots f 2 [D ' d \ D = 2,3,... that converges towards the actual value of u'(0)/2 
for a given value of d. We compare sequences with different values of d for 
inner consistency (all of them should give the same limit). 

Present approach is simple and straightforward: we just obtain the Taylor 
coefficients fj from the differential equations (J2J) or ^ in terms of f 2 , construct 
the Hankel determinant, and calculate its roots. Since f± is the first nonzero 
coefficient that depends on f 2 we choose Hankel sequences with d > 3. 

The Hankel determinant Hp exhibits many roots and their number increases 
with D. If we compare the roots of Hfy with those of Hf ) _ 1 we easily identify 
the sequence f? A that converges towards u'(0)/2. 

Present HPM may be considered to be a systematic generalization of Tu's 
approach [12] and is clearly different from the strategy proposed by Epele et 
al [13]. We stress that we have been using the Hankel condicion (fTUj) for quite 
a long time [22-30]. 

We first consider the TF equation for atoms. The convergence of the HPM 
has already been discussed in our earlier paper [21]; therefore, here we simply 
compare our results with those of Parand and Shahini [20]. Our estimate of 
u'(0) —1.588071022611375313 is by far more accurate than the one obtained 
by Parand and Shahini [20] -1.5880702966, Liao [18] -1.58606 and Khan and 
Xu [19] —1.586494973. Notice that all those authors kept several inaccurate 
digits in their results. This practice is misleading because it suggests that the 
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calculations are more accurate than what they really are. 



Table [T] shows our earlier results [21], those of Parand and Shahini [20] and the 
numerical calculation of Kobayashi et al [9]. Those HPM results obtained from 
a relatively small [5/8] Pade approximant are more accurate than the numer- 
ical ones for x < 1. We appreciate that the simple and straightforward Pade 
approximants are far more accurate than the more elaborated pseudospec- 
tral approach [20]. Notice that also in this case Parand and Shahini [20] kept 
several misleading inaccurate digits. 

Fig. ffl shows L D4 = log |2/j D ' d] - 2f [ 2 D ~ 1 ' d] \ for D = 3,4, . . . that provides a 
reasonable indication of the convergence of the sequence of roots for the TF 
equations for atoms in strong magnetic fields. We clearly appreciate the great 
convergence rate of the sequences with d = 4 and d = 5 (the rate of conver- 
gence for the case d = 3 is slightly smaller but still suitable for practical appli- 
cations). From those sequences we estimate u'(0) = —0.93896688764395889306 
that is exact to the last digit and considerably more accurate than the results 
published earlier [1,4,5]. 



4 Conclusions 



In this letter we have shown that the well-known simple and straightforward 
Pade approximants provide much more accurate results for the TF equation 
than the more elaborated Chebyshev pseudospectral method [20], exactly as 
we did earlier with the rather cumbersome HAM [21]. We also applied the 
HPM to the TF equation for an atom in a strong magnetic field, analyzed the 
convergence of the Hankel sequences towards the slope at origin, and obtained 
its value with unprecedented accuracy. 

The TF equation is an example of two-point boundary value problems that 
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are most important in theoretical physics. We have already applied the HPM 
to other such problems [31]. There are alternative approaches for the accurate 
treatment of two-point boundary value problems. Here we mention the work 
of Boisseau et al [32], Bervillier et al [33,34], and in particular a recent com- 
prehensive discussion of power-series methods for ordinary differential equa- 
tions [35]. 
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Table 1 

Values of the Thomas-Fermi function u(x) obtained by present method, the Chcby- 
shev pseudospectral one [20], and numerical integration [9] 



X 


HPM 


Chebyshev 


Numerical 


1 
1 


n /19/inns 

U.4fcZ4fcUUo 




n /19/im 


5 


0.078808 


0.078277758 


0.078808 


10 


0.024315 


0.025044744 


0.024314 


20 


0.005786 


0.006585633 


0.0057849 


50 


0.000633 


0.000761317 


0.00063226 


100 


0.0001005 


0.000023409 


0.00010024 




Fig. 1. L D4 = log 2ff' d] - 2fP~ lA 



for d = 4 (circles) and d = 5 (filled circles) 
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